clc; clear; close all;

%% Load data from the microstructure processing routines
load('intermediate\MicrostructureEstimatesAndOutputVaryingBetas.mat');

%Run file to input cross-currency basis data
run functions\CIPGraph3MEURCHF;

%Extract from the CIPGraph3MEURCHF processing the various cross-currency
%bases
dateVecCIP = CIPAllS1{:,1};
CIPCalcsAll = [CIPCalcAUD CIPCalcGBP CIPCalcCAD CIPCalcEUR CIPCalcJPY CIPCalcEURCHF];

%Plot the cross-currency bases
plot(dateVecCIP,CIPCalcsAll)
dateaxis('x',2)

CIPVecList = {CIPVecAUD,CIPVecGBP,CIPVecCAD,CIPVecEUR,CIPVecJPY,CIPVecEURCHF};
counter = 1;
dateVecCIP = dateVecCIP(3077:end,:);
dateVecCIP = dateVecCIP(selectorNonHolidays,:);

%Initialize matrix to store welfare and arb size calculations, matrix
%alternates between welfare and arb size computations, moving through each
%market (column 1 = AUD/USD CIP Welfare, column 2 = AUD/USD CIP Arb
%size...)
BetaValues = (0.3:0.025:0.6)';
WelfareArbs = zeros(305,12*13);
FunctionEvals = zeros(305,6*13);

%Initialize vectors with contract parameters given by contract
%specifications, including tick size (minimum fluctuation) and size of
%contracts (how many dollars a contract is worth typically during our estimation period)
TickSizes = [0.00005 0.0001 0.00005 0.00005 0.0000005 0.0001];
ContractSizes = [75000 81250 75000 150000 120000 150000];

%Compile coefficients for dailies
betacounter = 0;
while(betacounter<size(BetaValues,1))
    %Cycle through each contract and compute the welfare gains to arbitrage and
    %gap-closing size
    betacounter = betacounter + 1;
    time = 0;
    counter=1;
    while(counter<=5)
        %Extract the interest rates, spot, and forward data
        CIPVec = CIPVecList{1,counter};
        %Isolate only the recent period for which we have microstructure
        %estimates by day
        CIPVec = CIPVec(3077:end,:);
        %Create a selector vector that eliminates holidays and select out data
        selectorCIP = zeros(315,1)+1;
        selectorCIP([8 9 13 65 85 130 195 260 270 275],:) = 0;
        selectorCIP = (selectorCIP>0);
        CIPVec = CIPVec(selectorCIP,:);
        freeBetas = coeffs(:,(counter-1)*size(BetaValues,1)+betacounter);

        time = 0;
        %Cycle through each day in the recent microstructure dataset and
        %compute
        while(time<size(CIPVec,1))
            time = time + 1;
            %One can try varying the parameter to calculate welfare since this
            %is the cross-price impact assumption: robust=10 implies the price
            %impact is only 1/10 of what it normally would be due to the fact
            %that spot and forward markets feature significant spillove price
            %impact (going short one leg drives the other leg down and vice
            %versa, resulting in largely cancelling effects from the trade)
            robust = 10;
            robust2 = 100;
            %Calculate assuming (zero coefficient) no interest rate price
            %impact for simplicity since it is second order effect
            Values = @(m) [(CIPVec(time,1)+sign(m)*(TickSizes(1,counter)*freeBetas(time,1)/robust).*(abs(m)*1000000/(ContractSizes(1,counter))).^BetaValues(betacounter,1))./(CIPVec(time,2)-sign(m)*(TickSizes(1,counter)*freeBetas(time,1)/robust).*(abs(m)*1000000/ContractSizes(1,counter)).^BetaValues(betacounter,1)) .* (1./(1+(CIPVec(time,4)/100 - sign(m)*(0.0000000/robust2).*(abs(m)*1000000/120000).^0.5733))).^0.25 - (1./(1+(CIPVec(time,3)/100+sign(m)*(0.0000000/robust2).*(abs(m)*1000000/100000).^0.5852))).^0.25]';
            %Values2 is not used
            Values2 = @(m) [(CIPVec(time,1)+sign(m)*(TickSizes(1,counter)*freeBetas(time,1)/robust).*(abs(m)*1000000/(ContractSizes(1,counter))).^BetaValues(betacounter,1)) (CIPVec(time,2)-sign(m)*(TickSizes(1,counter)*freeBetas(time,1)/robust).*(abs(m)*1000000/ContractSizes(1,counter)).^BetaValues(betacounter,1))  (CIPVec(time,4)/100 - sign(m)*(0.000000013/robust2).*(abs(m)*1000000/120000).^0.5733) (CIPVec(time,3)/100+sign(m)*(0.000000/robust2).*(abs(m)*1000000/100000).^0.5852)]';
            %Set optimization options
            options = optimoptions('fsolve','FunctionTolerance',1e-12,'OptimalityTolerance',1e-24,'StepTolerance',1e-12,'MaxFunctionEvaluations',1e3);
            %Use solver to find gap-closing size
            [m1,m1Fn] = fsolve(Values,1,options);
            %Numerical integration routine
            grid = (1:((abs(m1)-1)/200):abs(m1))*sign(m1);
            evals = Values(grid);
            TotalWelfare = sum(evals)*((abs(m1)-1)/200);
            %Store computations
            WelfareArbs(time,2*((counter-1)*size(BetaValues,1)+betacounter)-1)= TotalWelfare;
            WelfareArbs(time,2*((counter-1)*size(BetaValues,1)+betacounter))= m1;
            FunctionEvals(time,(counter-1)*size(BetaValues,1)+betacounter)= m1Fn;
        end
        counter = counter + 1;
    end

end
mean3Ms = mean(CIPCalcsAll);
std3Ms = std(CIPCalcsAll);
min3Ms = min(CIPCalcsAll);
max3Ms = max(CIPCalcsAll);
meanabs3Ms = mean(abs(CIPCalcsAll));

%These are the data used for the CIP Calculations:
%ADBB3M CMPN Curncy Australian Bank Bill 3 Month
%BP0003M Index STERLING ICE LIBOR
%CDOR03 Index 3 Month Canadian Banker’s Acceptance Rate
%GTDEM3M Corp German Treasury Bill (Bubill)
%GB03 Govt (Three month T-Bills)
%JY0003M Index Japan 3 Month LIBOR


data2 = [mean3Ms;std3Ms;min3Ms;max3Ms;meanabs3Ms];
input.data = data2;
input.tablePositioning = 'h';
input.tableColLabels = {'AUD/USD','GBP/USD','CAD/USD','EUR/USD','JPY/USD','EUR/CHF'};
input.tableRowLabels = {'Mean 3M Cross Currency Basis','Std 3M Cross Currency Basis','Min 3M Cross Currency Basis','Max 3M Cross Currency Basis','Mean Abs 3M Cross Currency Basis'};
input.dataFormatMode = 'column';
input.dataFormat = {'%.4f',5};

save('intermediate\CompileAllCIPCalcsVaryingBeta');
